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Cosmological departures from general relativity offer a possible explanation for the cosmic ac- 
celeration. To linear order, these departures (quantified by the model-independent parameter rzr, 
referred to as a 'gravitational slip') amplify or suppress the growth of structure in the universe 
relative to what we would expect to see from a general relativistic universe lately dominated by a 
cosmological constant. As structures collapse and become more dense, linear perturbation theory 
is an inadequate descriptor of their behavior, and one must extend calculations to non-linear order. 
' If the effects of gravitational slip extend to these higher orders, we might expect to see a signature 

, of i:i7 in the bispectrum of galaxies distributed on the sky. We solve the equations of motion for 

■ non- linear perturbations in the presence of gravitational slip and find that, while there is an effect 

on the bispectrum, it is too weak to be detected with present galaxy surveys. We also develop a 
\ formalism for incorporating scale dependence into our description of gravitational slip. 
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I. INTRODUCTION 



A universe in which gravity obeys the laws of general relativity and is filled with baryons and cold dark matter 
ought to decelerate. If we write the spatially-flat Robertson Walker metric 

2 [ deceleration means that a = da/dr decreases with time. This is not the case in the Universe in which we live. 
' Observations of distant supernovae indicate that our universe is, in fact, accelerating (a is growing with time) [J, [l]- 
] To date, there are many theories vying to explain this acceleration. These theories can be generally divided into two 
categories: theories of dark energy and alternative theories of gravity. If we write Einstein's equations as 
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' dark energy attempts to explain the acceleration by modifying the right hand side (i.e., by positing that the universe 
, is filled with an exotic new material). Alternative theories of gravity attempt to explain the acceleration by modifying 
• ' the left hand side (i.e., by supposing that the laws of gravity obey different equations of motion on cosmological scales 
I than they would under general relativity). Great diversity exists within both categories. 

' Dark energy can be Einstein's cosmological constant |3|, a uniform energy density associated with the vacuum. 
[ Although this is the explanation currently favored by the data, there is no theoretical calculation justifying the 
observed value of the constant. Dark energy may also be a cosmological scalar field [3| evolving slowly through its 
potential V such that the equation of state w = p/p < a,t late times. If this were the case, we would expect to see 
w evolve with redshift at early times, which we do not. Dark energy could also be a vector field Q, though these 
^ ■ theories often introduce preferred reference frames or couplings between scalar and vector perturbation modes which 
we have yet to observe. 

Alternative gravity can be a scalar-tensor theory Q in which the dark energy scalar field is non-minimally coupled 
to the curvature terms in the Einstein-Hilbert action. These theories tend to predict departures from Newton's law of 
gravitation, though the simplest of them merely manifest themselves as a rescaling of Newton's constant G. There are 
also theories which modify gravity by introducing an arbitrary function of the Ricci scalar f{R) into the gravitational 
Lagrangian 0, H Q • These theories tend to introduce new scale-dependent effects into the evolution of large scale 
structure. Tensor- vector-scalar (TeVeS) theory [13, [HI adds tensor and vector fields to the mix of scalar-tensor theory. 
These theories require precise couplings to avoid the pit-falls of vector dark ener gy noted above. Multi-dimensional 
braneworld theories (like the Dvali-Gabadadze-Porrati - DGP - model) [T^. [l3l.Tl^ attempt to account for cosmic 
acceleration by allowing gravity to act in dimensions outside of our 3 + 1 Universe, imposing a scale beyond which the 
expected gravitational attraction is damped. Theories inspired by quantum mechanics, such as a Lorentz-violating 
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massive graviton [15[ , reproduce solar system tests of gravity but introduce new classes of cosmological perturbations 
which may or may not show up in our measurements of large scale structure, depending on initial conditions. 

All of these theories have solutions that can provide late-time cosmic acceleration. Successfully navigating this 
labyrinth of viable theories requires a set of observations complementary to the expansion of the universe for which 
different theories of gravity or dark energy predict different effects. The growth of cosmic structure provides one such 
set of observations. We model that growth using cosmic perturbation theory. If we consider only perturbations that 
transform like scalars, the first-order perturbed form of equation ([T|) is 

ds^ = a^[-(l + 2V')dT2 + (1 - 2(t))dx^]. (2) 

Neglecting the cosmic scale factor a, the potentials 4' ^nd </> - known respectively as the Newtonian and longitudinal 
potentials - supply the right hand sides of Newton's law of gravitation x = —Vip and the Poisson equation AnGa^Sp = 
V^0. Consistent with Newtonian dynamics, general relativity in the presence of non-relativistic stress-energy predicts 
that ^ = ^J. Alternative theories of gravity make no such guarantee. Each of the gravity theories cited above implies 
its own unique relationship between and -0 @, d, [^, [ll|, 0, 0, EB] • 

Lacking theoretical justification to prefer one alternative gravity theory over another, we will content ourselves 
searching for evidence of gravitational slip (our word for the difference between cf) and ijj) in general. Following 
Caldwell et al. p^ . we parametrize the relationship between cj) and "0 as 

= (1 + ^)0 (3) 
n7 = wqu'^ . (4) 

Obviously, the scale independence and redshift dependence of equation ^ are assumptions on our part. Because we 
arc interested in gravity theories which might explain the late time cosmic acceleration, we suppose that « ■0 at 
high redshift and that gravitational slip should grow as the inverse of the matter density. For the purposes of our 
calculations, wc will assume that the expansion history of the Universe exactly matches that of a ACDM universe, 
that is a universe which obeys general relativity (vo—O) and in which the acceleration is caused by a cosmological 
constant contributing a fraction J^a < 1 to the current density of the Universe. 

Equation ([3|) is obviously not the only possible parametrization of modified gravity. References [13, [H, [H, HO, HH, 



|22| . l23l l24j all explore different model-independent expressions of ^ 0. Reference [25| discusses the consistency of 
our choice Q with these other paramctrizations. In the end, it is simply a question of nomenclature. 

References psi . [26j , derived the linear-order equations of motion by perturbing a Robertson Walker metric ([T|) in 
the presence of a homogeneous, isotropic stress energy tensor 

= diag(-p,p,p,p) (5) 

and constrained wo against CMB, supernova, and weak lensing (of galaxies) data sets. Specifically, reference [26l | 
found the constraint wq = 0.09to;5g (2cr). Reference [l^l also tested equation ([3]), but against CMB and weak lensing 
(of the CMB) data sets. They also promoted the redshift dependence of equation ([4|) to a new free parameter. They 
found wq ~ 1.67l^ g7 (^c)- Both results are consistent with a ACDM universe obeying the laws of general relativity, 
but with significant room for departure (0 = "0 within a factor of a few). In this paper, we will attempt to complement 
those explorations by calculating the effect of non-zero tu on the growth of structure in the Universe beyond linear 
order. 

As gravitational collapse of structures progresses, the gravitational fields within overdense regions of the universe 
grow so that effects beyond linear order become important. Even in a ccj = universe, this causes the distribution 
of overdensities to evolve away from their initial Gaussian distribution p8| . These departures from Gaussianity are 
evidenced in the bispectrum, the Fourier transform of the three-point correlation function of galaxies. By altering 
the growth of structure [l^ , w ought also to alter these departures from gaussianity, an effect which we hope will 
be detectable in modern galaxy surveys. Section [IT] will derive the equations of motion for cosmic overdensities in 
the case of 7^ to "quasi-linear" order. Section IIIII will translate these results into Fourier space. Section IIVI 
will calculate the resulting effects on the bispectrum. Section IIV Al includes the effects of galaxy bias and redshift 
distortions. Section IVl will consider the effect of adding scale dependence to tu. Section IVll will discuss our results. 
Appendix |^ compares the results of our equations of motion (derived in Section |TT] using Eulerian Perturbation 
Theory) to previous results derived in Lagrangian Perturbation Theory. Appendix [B] discusses the averaged second 
order moment of the overdensity distribution (the skewness). Appendix [Cl extends our results to the next higher 
quasi-linear order, calculating the effect of scale-independent zu on the kurtosis of the overdensity distribution. 
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II. EQUATIONS OF MOTION 



This paper will work in the convention of Eulcrian Perturbation Theory (see Appendix |^ for a discussion of 
Lagrangian Perturbation Theory). We will directly evolve the perturbed quantities (p, 5 = {p — p)/p and zT, the 
perturbed velocity field in the cold dark matter fluid. References [igI . [25l [26j discuss how to modify the linear-order 
equations of motion for these perturbed quantities in the case oi w ^ Q. For the sake of brevity, we will merely 
summarize that discussion here. Of the linear-order perturbed Einstein equations presented in reference p9| . the 
time-time and diagonal space-space are discarded. The time-space Einstein equation is preserved as a consequence of 
the assumption that the cold dark matter fluid remains, on average, at rest in our reference frame. A term is added 
to the off-diagonal space-space Einstein equation so to provide for n7 as a new source of cosmic shear. Because we 
have discarded the time-time and diagonal space-space Einstein equations, we can no longer assume that the Poisson 
equation is valid. This will be important when deriving the quasi-linear equations of motion below. 

At quasi-linear order, gravitational collapse has advanced such that 

,5,w/c> 0/c^v/c^ (6) 

i.e. the dynamics of the local fluid dominate the dynamics of the spacetime. The metric is still Robertson Walker, 
but now, the stress energy tensor ([5]) is replaced by that of a perfect fluid 

T^"" = {p + PW'u'' + g^"'p 

where u^^ are the components of the fluid's four- velocity [30| . To consistently determine the smallness of perturbed 
quantities, we explicitly include factors of c in the perturbed metric giving 

We derive our equations of motion by solving for the dynamics of this reformulated system to zeroth order in 1/c. 
This is equivalent to assumption ([6|). 



Requiring u^u^ = — 1, we find that 
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where 7 = 1/-^/! — iP- jc^ , as usual. Thus, to the required order (because we are considering perturbations in the 
matter distribution, we set p = 0), 



V^,T^^° ^ ^{5 + d,{l + 6y)=Q (9) 



V^T^* = ^ {v' + Hv' + v^djv' + diil>) ^ (10) 



a2 

= ~8TrGa^v'p(l + S)^ (11) 

where H, is the conformal time Hubble parameter 71 = a/ a. Note that, to lowest order, 9' = a^^S^^dj. Equations ([9|) 
and (fTO|) also correspond to equations (2) and (3) of Catelan and Moscardini's paper deriving the quasi-linear fluid 
equations of motion in unmodified general relativity [sij . Their third equation of motion derives from the Poisson 
equation, which we discard in favor of our equation (jlip . (Note that Catelan and Moscardini use the coordinate time 
t where we use the conformal time r). Following Catelan and Moscardini's lead, we substitute equation ^ into the 
divergence of equation (fTO|) to get 

5 + H5 = {1 + zu)diSdi(l3 (12) 

+ay [(1 + syv^] + (1 + ,5)(i + 

where we have used equation ([3]) to rewrite "0 in terms of (p and za. In our notation = didi. Using equation ([9]), 
we can rewrite equation (jlip as 

^n^n^6 = (j) + H{1 + w)(l)^ (13) 
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the time derivative of which is 

3 
2 



Using equation (|T3|) . we find 
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^n,nn^(d + HS^ = y^i^4> + H{3 + w)^ (i4) 

Equations (jl2p and (|14p provide an algorithm by which we can solve for (p and 5 to arbitrary order. Once we have 
solved for and 5, it is a simple matter to use equation ([9]) to find v. 

Assume that (/>, 5 and v can be expanded as </> = etc., where (/)^*^ » (i.e., "0^*^ is the ith order part 

of (/)"). In that case, we can expand equations (fT2|) and ((T4|) to a given order n, then use the proportionality between 
their left hand sides to get a single equation of the form 

^0^") + <^'"'H(3 + w) + [^'3(1 + - Qra) + Htu] j = ^r!,„H25(") (15) 

5^") = ^ |(l + tz7)5,<5('^)a,<^(''' +5„-[(l + (5('^))u('')''w(^)^'] +(5('^'(1 + ^)V'0^'')} (16) 

a+b+c— n 

where the source terms 5^"^ come from the non-linear part of equation (|12p . After solving equation psp for V^0'"\ 
one can use equation ([T! 



^5'(«) + (l + ro)V2</,(") (17) 

to solve for (5^"-'. 

Comparing equations ([T5|) and P7|) . we see that, even though the Poisson equation no longer holds, (j) and 6 are 
still separable at first order (0 = f{T)(p{x) and 5 = D{t)^{x)) and related such that the spatial parts of S are the 
Laplacians of the spatial parts of (i.e., ^ ~ V^i^s). Indeed, to first order 

= /(r)v.(f) (18) 
= L'(t)vV(:?) (19) 
u'l) = -D^ip (20) 

/ + /W(l + tu) = ^a,nH'-D 

b + nb ^ {i + w)f. 

At higher orders, the expressions are less compact 

0(2) ^ a(T)^(f ) + /3(r)S(x) 

^(2) ^ v^V^A + V^V^B (21) 
1/(2) = VA(-i)" +£>£))- VSX)^ (22) 



where 



V^A a, (VV9,;(^) (23) 

V^B = a,, {d^^djip) (24) 

a = -d7i:(3 + vj)+ '^^^'H^{\ + n7)L>/ - a [m^{\ + - r2,„) + wU] (25) 

/3 = -/37^(3 + n7) + |l]™W2£)2_^^3^2(^^^^^^_^^-)_^^-^j ^26) 

= -nif + {l + vj)a+{l + vj)Df (27) 

i)^ = -nv'' + {l + m)P + b^. (28) 
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In the w = limit, these expressions agree with the GR results of [S^- There are seven different spatial functions 
comprising (f)^^^ . We present them in Appendix IC II 

We now only lack initial conditions on the relevant growth terms in our attempt to solve for 6 in the quasi-linear 
regime. Since the power spectrum Pg is initially Gaussian, non-linear growth terms (i.e., a and f3 from equation (|2ip ) 
are integrated from the early-time initial conditions = 0, di = 0. To find the initial conditions for the linear growth 
terms / and D, we combine the first order equations of motion (|f 2|) and (|13|) . giving 



D + DH = {l+rjj)f 
3, 



f + fHil + w) 



to get 



-n„,H'D -H'D -no 



I 



-HD - HD 
d 



-HD + C 



(29) 



where C is a constant. Since we assume that lima-,o n7 = 0, we use the ACDM equations of motion combined with 
the Poisson-equation result / = (3/2)OmW^-D to find the correct value for C for a given fi^. This value changes 
depending on the initial conditions assumed for D and D. However, the normalization of the measured statistic (j56p 
means that our final results are resilient to such choices. Figure [T] plots the effect of gravitational slip on the first 
order growth function D and its derivative. Figure |5| plots the same effects on the second order growth functions 
and V''. 
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FIG. 1: We plot the change in the first order growth functions D and /i = resulting from varying zuq. Note that /i 

is a different parameter from the / used in equation (|18p . All other parameters are set to be the WMAP 5-year maximum 
likelihood values [U- The change is calculated relative to a ci7 = (GR) model so that AD/Dgr — {D ~ Dgr) / Dgr- As 
found in [2^, -ujq > amplifies the growth of 5^^\ while vjq < suppresses it. 



III. FOURIER TRANSFORM OF {S") 

In Section ITVl we will derive an expression for the bispectrum of the matter distribution in the case of 7^ 0. This 
calculation will require taking the Fourier transform of terms like (i5"(a;)) for n > 1, which we illustrate below. 
If we expand 6{x) ~ J2i then 

= ((<5(l))")+n(j(2)(5(l))'-l)+„(5(3)(^(l))n-l)+ Q((j(2))2(^(l))n-2^^ (gg) 

(note: superscripts in parentheses are orders of expansion; superscripts outside parentheses are exponents). Therefore, 
in order to get the Fourier tranform of ((5"), we first need the Fourier transform of ((S'*-*). From the form of equation 
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FIG. 2: We plot the change in the second order growth functions I>" and D'' resulting from varying zuq- All other parameters 
are set to be the WMAP 5-year maximum likelihood values [sj- In the case of solid lines, AC = — Cq^. In the case of 
dashed lines, the second order growth functions are normalized by a factor of D'^ {D is the first order growth function) so that 
AT?* = V (Dor/ D)^ — 'D'qr- The normalization greatly diminishes the effect of varying zuq, indicating that most of the effect 
of zajq on the higher order growth functions enters as a normalization. 



(fT7)) . we can see that 

<5«(^) = i?«(r)A(v.i,^2,...,^0 (31) 

where Z?'*-* stands for some combination of ith order growth functions and A is some combination of spatial derivatives 
acting on the i powers of ip{x) (see, for example, equation ET|) . Taking the Fourier transform of this qualitative 6^"^^ 

we find (here we will introduce the notation Sa'' = 5^'^^{ka)) 

^^^^(^) = / (^^e-''^^WA(^i,...,(^0 (32) 



(PxSki 



_i_le-(S5..'^.-^)^W/C[]^,„ (33) 



(27r) 

-A|_^,3^(fc_^,,)i^«/Cn^™ (34) 

(note d^ki i = 11}=! d'^kj) where, in going from equation ([5^ to ([33]), we have taken the Fourier transform of each 

individual factor of ip{x) in equation pip . /C represents the combination of wave vectors {fci, fc2, . . . fc^} deriving from 
the differential operator A. 6^ is a Dirac delta function. We used equation (|19p to go from equation ([31)1 to ([35]) . 
From equation ([35]) . we can see that a term of the form (5^*^(.t)((5^^^ (2?))" can be written in terms of Fourier transforms 

as 



<5«(.f)(<5(l)(f))" = / f^f!^^^Sik+J:-^,k,) ^ TT -gil) 

V"^ / 771—1 



^ j d^k<fk,...n d^k[ .s(Uy.U fe^.) TT ja) ^5Uk-y k[) ^7^!''^^." ^ fr C (36) 
J (27r)i("+i) (2^)i(»-i) " ^ {D^^^Y f;}^ k'i^ ^ > 

where, to be explicit, the k'- wave vectors come from the Fourier transform 5^^\ Integrating the right hand side of 
equation ([36]) over the wave vector fc, we see that 

jW(f)(5W)" = /^!^l^e"-(^S"^"^)x^^— ^ff^W 
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whence 

(<5«(f)(<5(i)r) = / i^%±!l-e"-ES?^^) X ^7^^'^^'^ ^ X (^\..^^. (37) 
\ ^ A ; / y (27r)i^+") {D^^^Y HLi 

At this point, it is useful to recall that, for a Gaussian field (like (J*-^-*), 

{5^^^5'iK..~5\l^) =0 nodd (38) 

+ permutations n even (39) 

so that, for example 

= 

From this, we see that the right hand side of equation ([57)1 involves {i + n — l)(i + n — 3) . . . (1) terms if {i + n) is 
even. If {i + n) is odd, equation ((37)) is identically zero. Finally, we recall that 

{5'^^^~5^^)=5l{k,+k2)Ps{ki) (40) 

Combining these results with equation (|30p gives us an algorithm for evaluating the higher order correlation functions 
of the 5 distribution. We illustrate this by considering the bispectrum. 



IV. THE BISPECTRUM 



The bispectrum is the three point equivalent of the power spectrum, i.e. it is the Fourier transform of the three-point 
correlation function. In the expression 



{5{x,)5{x2)5{x^)) 



(2^)9/2 



1=1 



(27r)9/' 



3 »E;U^,fc.S(fci,fc2,4) 



B is the bispectrum. From equations (l30|) and ([38]) . we see that the three-point correlation function will, to leading 
order, be made up of three terms like 

(5(2)^(1)^(1)). 

We write S^'^\ using equations ([221) and (HH), as 

j(2) ^ V''(^^^(p^iV^(p + V^ipW^^'^ +v''(^i■i(p^^j(p + \/'^^V'^(p + 2^iip^iV'^ip 
From equation (|35p. we see that this means 



|(2)(fc) 



d^k'd'k"^-^WLlsl{k -k'- k") X 



2?'^ 1 



k' ■ k'' 

fc'2 



1 + 2 



k' ■ k" {k' ■ k"Y 

fc'2 ^ k'^k"^ 

'/\\2 



d^k' 5{k')5{k-k')T / k'-(k-k')\ k'-(k-k') fc'- fc-/c'M , ^ 

(27r)3/2 1)2 L V A:'2 y V fc'2 y2(i„_y^2 )\ ^ ^ 



k'-^{k-k'Y 
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We can now average equation (|4T|) with two factors of S{k) to find the form of the bispectrum. 



D 



^4 



2?" 1 + 2 



k4 ■ (ks - ki) [ki ■ (ks - ^4)) 



2 V n 



^4 



kiik^-k^y 



D2 



V 1 



■{P(fcl)P(fc2)<5|5(i+3 -4)'^i)(2+4) +4)^i>(2+3-4) 



A;4 • (fcs - ^4) 



^4 



p" 1 + 2 



fc4 • (fcs - fc4) (fc4 • (fca - ^4))^ 



A:|(/c3-fc4)2 



(42) 



(note ^|)(j^j„;-, = 6^{ki + kj — fc;)) where we have used equation (|40p to get the power spectra and factors. 
Note that, by the coefficients of the growth terms, the S^{k3) term will vanish upon integration. Once that term is 
discarded, we see that the bispectrum is proportional to S^{ki + k-i + kj,). Specifically 



B{ki,k2,k3) = 



P(ki)P{k2) 



£12 

I P{k,)P{k2) 
+permutations 



^2 • ki 

1.2 
^2 



)+P''(l + 2 



fc2 • fcl ik2 ■ fcl)^ 



k^kf 



6Uki+k2 + k3) 



P°(l 



k2 ■ ki 



)+V\l + 2 



k2- ki (fc2 • kiY 



kl 



1^2 '^1 



5%{ki + k2 + k3) 



(43) 



It is straightforward to integrate the equations in Section |TI] and determine how zdq affects equation (|43p . It is less 
straightforward to turn this calculation into a real- world constraint on voq. 



A. Galaxy bias and redshift distortions 



It is presently impossible to measure the dark matter density at all points in space. As its name suggests, we 
cannot see dark matter. We can only see the galaxies that form in dark matter haloes. Going from astronomical 
observations to a determination of the bispectrum (j43[) requires assumptions about how the galaxy distribution tracks 
the dark matter distribution (the "galaxy bias" ) for which we have little theoretical motivation. As if that were not 
hard enough, we only actually see the galaxies in two dimensions (altitude and azimuth relative to our telescope). 
We infer the radial distance to galaxies by measuring their redshift and assuming that Bubble's law is valid. This is 
a decent assumption for redshifts of a few. Unfortunately, it means that galactic peculiar motions interfere with our 
determination of galaxies' positions, giving rise to "redshift distortions" in the observed distribution of galaxies. It is 
known how to correct for these effects in calculating equation We will do so below. Our derivation relies heavily 
on the tu = calculations presented by Bernardeau in Section 7 of reference [s^l • 

The correction for scale- and time-independent galaxy bias is simple. Assume that the excess number density of 
galaxies 5g — {ug — fig) /fig relates to the matter overdensity 5 by 

i 

where bi are (constant) coefficients of the expansion. If we then write the overdensity S as we did in the discussion 
leading up to equation ([30|l . we have, to second order 

53 = 6i<5W+6i<5(2) + |(,5«)2 + ... (44) 

To first order, the expansion is simply dg = biS. 

To account for redshift space distortions, we follow Section 7 of reference [s^] or Section 9.4 of reference [1^. We 
denote redshift space (in which radial distance is reckoned from Bubble's law) by Xs- Physical space will remain x. 
Because observations are made in redshift space, we want to calculate the bispectrum of the redshift space distribution 
using the physical space evolution equations we derived in Section |TT1 We will work in the plane-parallel approximation 
in which the sky is flat and the radial direction is z. By Bubble's law, the redshift space z coordinate of a given 
galaxy will be 



(45) 
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We use the conformal time Hubble parameter because, throughout this work, v has also been calculated as the 
conformal time velocity. It is also true that Tip = ^^o s-s long as Cq = 1. From equation (|45p . we have 



dxgdysdzs -, 1 ^ 

1 + Trr"^^^- 
rto 



d^x dxdydz 

Because we do not want our change in coordinates to change the mass of a given region, we set 

{l+6s)d^x, = {l + 5)d^x 



(46) 



from which we find 



d^x 

{l + S)^„ 1 



d^a 



1 + 6 



d'^Xs\ d'^x 



S dzVz 

no 



d^x / d^. 
dh 



d3a 



Now, we can write 



Ss{k) = 



(27r)3/2 

d^Xs 
(27r)3/2 

d^x 
(27r)3/2 



^'S,{Xs) 



Tio 



fx_ 

d^Xs 



— ik-x 



S - —dzVz 
rto 



(47) 
(48) 



where the « comes from expanding equation (|47p to first order in perturbed quantities. If we want the Fourier 
transform of the galaxy overdensity in redshift space, equation (|48p must be rewritten 



dh 



(27r)3/2 



1 



Ho 



2 tio J 



(49) 



To proceed further, we will require expressions for the first and second order parts of the peculiar velocity field, if-^'^ 
and w^2)^ terms of the Fourier transform (p. 
Recalling equations ([20|l and ([22|) . we write 



^1) 



t(i) = 



-iD 



d^k 
(27r)3/^ 



e"^-^fc(^(fc) 



D k ~, 



t(2) 



^(2) = , 



{V" - DD){klk2) + V\[kl ■ k2]k2 + fc2fc2) 



= -^Dk^ = '^^s^'^ 

d^k' r ^^^^ 

■^'^{k - k') X {V - DD){k'{k - k'f) + V\k'[k' ■ {k - k')] + k'{k - k'f) 



(27r)3/2 



(50) 



(51) 



Hoping to find a redshift space expression equivalent to equation 
order parts. For the first order part, we find 



d^xd^k' 
(27r)3 



e"-('='-^)(-A:'2)D(^(fc') 



we separate equation (j49p into first and second 
d^xd^k' 



(27r)3 



„ix-(k'-k) {k'zY 



D(p{k') 



D 



= {-k%^D~—ki)m 

rto 



(52) 



where the factors of k^ come from dzVz (recall that Vz oc kz by equation [SO]) 



10 



The second order part gives 



D 
Ho 



62 
' 2 
(Pk' 
(27r)3/2 % 

d^k' 



d^k' 
(27r)3/2 



D^k'^{k~k'y<f{k')ip{k-k') 



(p{k')(p{k - k') X k'^{k - k')^{V'' - DD) + V\k' ■ {k - k')k'^ + {k- k'fk'^) 

b 



^^-^^{k')^{k - k') X [ - h.kXik - k'?D - Ji^^kXik. - k'^f 



(53) 



where 5^'^\k) on the right hand side of the first hne is given by equation (|4ip . Combining equations (|52|) and ((53)) we 
can now write, by analogy with equation 

f(l)j(l)x(2) 



d^ki{p{ki)P{k2)6%{ki + fcg - k4)6^D{k2 + ki) 

+P{ki)P{k2)S%{h + h)S%{k2 + ks- ki) + P{ki)P{ki)dUki + k2)SUk3)} 



2 £12^0 

b 



ki ■ {ks - ki) 



^4 



k^-ih-k^) {k4-{ks-k4)y 



kjiks-kiY 



biD- 



^{V^ Db) + P^(fc4-(fc3-fc4)fc4. 

k^ 

kszkiz 



kjiks-ki)^ 

kiz)\ 



kl 



no ki{k^-kiY 



hi 



D 



kiz 

DHo V kl 



hi 



D 

Imo 



(54) 



The galaxy bispectrum in rcdshift space Bg^s has the same form as equation (j43p with the appropriate substitution 
from equation Note that, though the b2S^{k^) term is not identically zero (as the (5|)(fc3) terms in previous 

expressions are), we are still justified in discarding it, as we will not be interested in values of Bg^s for which k^ = 0. 

The introduction of redshift distortions into equation breaks the isotropy of the B. Peculiar velocities only 
affect our measurement in the z direction, so it matters which way our {ki, k2, k^} triangles are oriented. We can 
average over this orientation-dependence by integrating 



Bi,vgiki,k2,k3,9i2) 



sm^i 
An 



dOi 



27r 



d(t>Bg^s{ki,k2,k2.) 



(55) 



where 9\ is the angle that ki makes with the z axis, (p is the angle that the plane defined by ^2 and ^3 makes with 
the plane defined by fci and the z axis, and 612 is the angle between fci and ^2. Appendix 1X1 checks the consistency 
of our expression with prior results derived in Lagrangian perturbation theory. 

To eliminate effects due to initial conditions, it is typical to talk about the normalized bispectrum Q defined so 
that 



Q(fci, fc2, k. 



6avg(fcl,fe,fc3,^'l2) 



(56) 



h\al[PiP2 + PiP3 + P2P3) 
P. = Ps{k,) 

2 1 o 

_ a dD 

" bli^' 

Figure [3] plots (Q — Qgr)/Qgr for different values of tuq and Vim- Qgk is the value of the normalized bispectrum ((5^ 
in a -07 = universe with the WMAP 5-year maximum likelihood cosmology [3^. We use a Monte Carlo integrator 
based on the sobol sequence generator sobseq() presented in reference [33| to evaluate integral ((55)) . The WMAP 
5-year team [s^ reports fi^ — 0.26 ± 0.03 (Icr). Reference [2^ finds that CMB data alone gives Wq = 1.712'° (2c)- 
We find that the variation in the bispectrum due to a Icr change in f2m is comparable to the variation due to a 
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FIG. 3: We plot the effect of varying u^o and n,„ on the unbiased normalized bispectrum with redshift distortions, Q, defined 
in equation (|56p . The vertical axis is the relative deviation from the result in the WMAP 5-year maximum likelihood GR 
cosmology. The horizontal axis is the angle between ki and k2 in units of n. Linear power spectra were calculated using the 
code CMBfast [i^ modified as in reference [2^. We find that Q is much more sensitive to small changes in Qrn than it is to 
large changes in zuo- This rules out the possibility of constraining -ouo from the three-point correlation function. 



2(7 change in vjq. Thus wc conclude that the three-point function will be of little help in constraining the value of 
TUo- Even if future experiments improve our constraint on the fact that measurements of the galaxy three-point 
correlation function are only precise to > 10% (see Tables A1-A3 of reference [s^) means that we are a long way off 
from being able to constrain -cdq from the distribution of galaxies in the Universe. 



SCALE DEPENDENCE 



As was discussed in reference [25| , our scale-independent model of gravitational slip ([4]) does not affect the shape of 
the power spectrum Pg. Any effect on the bispectrum observed in Figure[3]must therefore be due to the renormalization 
of second order perturbations relative to first order perturbations in the presence of w ^ 0. Clearly, this effect is 
weaker than that of actually changing the shape of Ps, as varying does. Figure [5] plots the effect of varying vjq on 
both the unnormalized second order growth functions and (solid lines) and the normalized second order growth 
functions P"/_D^ and V^/D^ (dashed lines). The normalization greatly reduces the effect of wq, lending credence to 
our hypothesis that the principal effect of scale-independent gravitational slip on the bispectrum is through the same 
renormalization previously found for the first order growth function (26j . Fortunately, the most popular alternative 
gravity models all predict scale-dependent effects [l9|. It is therefore incumbent upon us to consider the effect of 
including scale dependence in parametrization Q. In this section, we find that scale-dependent tu does amplify 
non-GR modifications to the bispectrum, provided one considers the right combination of {fci, fc2, fcs}. 



A. Scale-dependent equations of motion 

Because we want solar system tests to remain consistent with general relativity [s^ , we will work in Fourier space 
and impose the constraint that — > as fc ^ oo. In this case, it is an easy matter to calculate the effect of 
scale-dependent zu on our first order results, provided that we rewrite parametrizations ^ and @ as 

'4> = {l + w)(j> (57) 
zu = woK{k)a^ (58) 

with the restriction that lim^^oo K{k) = 0. 

To first order, the equations of motion for {(f>,S,v} remain unchanged under this reparametrization. The only 
difference in our results is that, under (|58p the growth functions D and / are dependent on the modulus of the 
wave vector k as well as on redshift. To higher orders, the presence of the non- linear source term (jl6p forces us to 
reformulate our approach towards calculating the bispectrum. 
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Our ability to write equations of motion (|251) - ([28|) depended on the separability of equation ((T5|) into scale- and 
redshift-dependent parts. Using the first order solutions for {0, S, v} we were able to write the second order source 
term (fTBl) as 



5(2) = D{1 + zu)/9,(^VV) + D^d,j{d,ipdj^). 
Under reparametrization (|57|) . we have 

+ b{ki)b{k2)[klkl + (fci ■ fca)' + • fc2 + klki ■ fcs] }(^(fci)(^(fc2). (59) 

The integration over wave vectors (stemming from the fact that reparametrization (j57p is defined in Fourier space) 
makes it impossible for us to find a separable solution for {(f)^'^\5^'^\'if^^}. This docs not, however, mean that we 
need be stymied in our attempt to calculate the bispectrum under (|57|) . 

Upon inspection of equation (|54p . we see that the bispectrum accounting for redshift space distortions and galaxy 
bias is made up of terms of the following types 

(60) 
(61) 
(62) 

(9,^(2)^(1)5(1)) (63) 
(9,i)(2),5(i)a,z)(i)} (64) 
(5,«(2)a,^(i)a,i)(i)) (65) 

(66) 

(5(i)<5(i)a,w(i)a,{)(i)) (67) 

(9,{}(i)a,v(i)9,«(i)9,{}(i)). (68) 

As noted above, we can calculate the terms built totally out of first order pieces (|66l67l68p exactly as we did in Section 
IIV Al allowing for the new scale dependence of the relevant growth functions. We can calculate the terms with second 
order parts by modifying the equations of motion and PT|) so as to calculate the ensemble averages (5(2)(5(^)(5(^)), 
etc. directly (rather than calculating the second order part, multiplying by the first order parts, and then taking the 
ensemble average). We explain this modification below. 

Consider equation ([TT]). If we multiply both sides by two extra factors of 5^^^ take the Fourier transform and take 
the ensemble average, we have 

d{Fouricr}|(^f)5('^4'^) = j d{Fourier}{(5p^5(i)4'^) - fc2(i ^ ^)(^(2)^(i)^(i)^| (gg^ 

where the integral over djFourier} stands for an integral over 

(27r)6 

which reduces to 

— i-3^(5|, (fc - ki - k2 - ks) (70) 

[ZTT ) 

upon integration over (Px. Simple differential calculus allows us to write the integrand of equation (j69p as a second 
order differential equation for 

(^(2)^(1)^(1)) 

which is exactly the form of terms like (j60p in equation (j54p . For example, 

isr'm'') = |:(^T^'5(i)4^))-(5T^^r^^^^V(^p^^^^'^^^^) 
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where we have used j'""^' = ~klD{ki)(p{ki) . We can use equation (fT5l) to write a similar differential equation, which 
we can evolve to calculate the (0(^)(5(^-'(5'^^') term on the right hand side of equation As for the source term in 

(|69)) . we can use equation (|59p to write 



4 

E 



(2^)15/2 "15V- - -jy^ ■ • • / " DM^kf ' ' V"^(l+3)"^(2+4) ^ "13(1+4) "15(2+3) 

where {. . . } represents the combination of k vectors and growth functions in the integrand of equation (|59p sym- 
metrized in terms of fci and fc2- The S^{k) in equation ([72]) means that the right hand side of equation (j69p will 
be zero unless fc = 0, which, from the form (|70|) of the integral measure in equation ((69|) . implies that the left hand 
side will be zero unless fci + ^2 + ^3 = 0. This is simply the familiar result that the bispectrum is defined only for 
triangular configurations of k vectors. Noting that 



{m''oAi')--^'hs?mi'') (73) 



we are now able to solve for the terms (|61l) and (p^ in equation 

To calculate the terms in (|54p which involve SzU^ , we use equation Q to find 

We can treat the first term on the right hand side of (|7i|) the same way that we treated terms in equation 

(I69l). Note that 



((v-)^f )^T^f ^) - - ^ + ^ 1 X ii^-rsrsi'W)- (75) 




The term ((V-i)(5f ^(5^ '4^^) 

can be found by using the same evolution equations we used to find ((5(2),5(i'j(i)} except 
eliminating the leading factor of 2 from the source term ([72]) and rewriting the {. . . } as 

{...} = D{k2)fik3){l + mik3))klk3 + Dik2)Dik3){klk3 + k2-k3k3). (76) 

The factor of 2 is eliminated because 6^'^ S^^"^) is not symmetric in its wave vector arguments. We will 

need to account for this by summing over all possible arrangements of {fci, fc2, ^3} when we calculate the bispectrum. 
The new scale dependence (j76p is to account for the inverse V operator. Once the equations of motion have been 
integrated, we must also multiply by a factor of —^12^32/^3 to account for the fact that we are interested in the dz 
derivative of the z-component of v^'^K 

The second term on the right hand side of (|74p is found (relatively) simply from the first order solutions for v and 
S. Again, we ultimately want to find a term that looks like {dzV^JS^^Sl^^^). Our first order solutions give 

™z{si%g^) ^ ^ I '^^SUk-k,-k,){k,z + k,M%g' 

- J (27r)3/2 «5jlK42 +K5zj-^-^04 "5 

where FT denotes a Fourier transform. From here we see that ¥i:{dz (<5(i)w(i)) (5(i)(5(i)) wiU contribute terms like 

□ D k3zkiz D3 
ki D3 
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to the bispcctrum. This contribution can be added to the contribution ^ to find that, the {d^v^S^^^S'^^^) term in 
the redshift distortion-corrected bispectrum will be 




(77) 



where the sum is over the 6 permutations of {/ci, fc2, ^3}. This is the term in equation ([5^ . Using equation ((75)) . 
we can also evaluate terms (|64p and (|55|) . Thus, we are able to calculate the redshift distortion-corrected bispectrum 
in the case of scale-dependent gravitational slip. 

Because our source term ((59)) depends on the power spectrum at specific length scales we cannot, as in Section HIl 
evolve our equations of motion (j69p from arbitrarily early times without regards to the shape of {4>^^\ S^^^} and then 
add the scale-dependence later. We must, instead, begin with realistic initial conditions for the perturbation modes. 
Since the form of ([55]) is such that lima^o nj = 0, we take these initial conditions from the Boltzmann code CMBfast 

po| with the appropriate background parameters and modified it to output {<f)^^\ (f)^^\S'-^\S^^'>} at a = 0.1 (chosen 
as an epoch early enough that tu and late enough that quasi- linear structure growth can be expected to set in). 

We still assume that 0'^^ = i^'^^^ = S'^^^ = = at this initial epoch. With these initial conditions, we can calculate 
the redshift distortion-corrected bispectrum ([5^ in the case of scale-dependent tu. All that remains is to explore the 
effect of parametrization ([57]) on obsevablc statisitics. 



B. Effects on observables 



To compare the effects of scale-dependent gravitational slip to a = universe, we perform the same average over 
orientation as in equation (|55p and use the same Q statistic as defined in equation (|56p . The denominator of Q in the 
case of scale-dependent gravitational slip is written 

b'i{PiP2aia2 + PiPsaia-s + ^2^30203) 

with 

1 + lfu^ + In 

a dD{ki) 
D{ki) da 

Figure m plots the effect of scale-dependent w on the Q statistic. For the purposes of this plot, we choose vu to have 
a scale dependence given by (recall equation [55)) 

K{k) = ^^y^ (78) 

In Figure [4] we choose fccrit = 0.01(Mpc)~^ and zjq = 5. This is an illustrative model only and is not meant to 
represent any specific theory of gravity. Including the scale dependence ([78]) amplifies the effect of gravitational slip 
on Q, but only for configurations {fci, fc2, fc3} in which the length scales correspond to the peak of K{k). It therefore 
seems reasonable to conclude that, while the bispectrum may be able to tell us something about scale-dependent 
gravitational slip, it will only do so if we choose to examine a particularly sensitive range of k. 




VI. CONCLUSION 



Contrary to expectations, wc have found that quasi-linear order structure growth is ineffective at constraining 
scale- independent gravitational slip ([4]). However, with judicious triangle choice, wc may be able to say something 
about a scale-dependent gravitational slip ([55)) by measuring the bispectrum. To date, there has been no attempt 
to constrain scale-dependent gravitational slip with even linear-order cosmological data. Conventional wisdom has 
been that the data is still too imprecise to say anything meaningful about so detailed an effect as scale-dependence. 
The next generation of cosmological experiments will hopefully remedy this oversight p6j . Given the prevalence of 
scale-dependent effects in alternative gravity theories [l9|, it may behoove us to find a useful, model-independent 
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FIG. 4: The same as Figure |3] except for scale-dependent vj of the form (|58|) with tuo = 5 and K{\z) given by equation (|78p . All 
length scales are in units of Mpc. All triangles considered are isosceles (i.e., fc2 = fci). As in Figure[3l Q is the angle between fci 
and ki. Figure [4 ( a) | plots K{k). Figure p:(b) | plots the change in the bispectrum relative to the WMAP5 maximum likelihood 
GR cosmology |36l| for different triangles. The case of fccrit = 10^'' is shown to illustrate what the curve would look like for 
scale-independent zu (|4]). The effect of rzr is much more pronounced for triangles whose sides correspond to length scales at 
which K{k) peaks. 



parametrization of the scale dependence K{k) (j58p and subject these future datasets to the same analysis given 
in reference [26l | . Depending on what such investigations say. the results of Section |V] may be able to extend our 
knowledge by bringing to bear measurements of the three-point correlation function on the question of gravitational 
slip. 
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APPENDIX A: COMPARISON TO LAGRANGIAN PERTURBATION THEORY 



Equation (I54J) was derived using Eulerian perturbation theory (i.e., we evolved the perturbed fields {(f>, 6, v} directly). 
Hivon et al. 41 1 perform the same calculation using the Lagrangian perturbation theory formalism laid out in 
references 4^ 43, H,!!^. Rather than evolving the perturbed density and velocity fields, Lagrangian perturbation 
theory evolves the displacement field ^'(t, q) defined as the displacement of a fluid element whose initial position in 
space was q. Put another way, at some time r, the position of a fluid element originating at spatial position q and 
initial time r, is 



X ^ q + *(t - Ti,q). 
From this, we find that the perturbed velocity field of Section [IT] is just 

« = — = *(t,(7) + («• V)«'(t,(7) 
dr 

= ^(1) + (^(1) . v)*(i) + ^(2) + 0(3) 



(Al) 



(A2) 



where we have set = 0. Note that in our notation (as in Section |TT]) , an overdot denotes partial differentiation with 
respect to conformal time r. In Hivon et aFs notation, an overdot denotes the Lagrangian derivative (dr + u • V). The 
same conservation of mass considerations that prompted us to write equation (j46p . now imply 

p{l + 6)d^x = pd^q. 

There is no S on the right hand side because we assume that the fluid starts from a homogeneous state. Algebra and 
equation (|Aip give 



fq 

d^X ' 



1= ( i + v-*(^^ 



(V.^(l))2 



(1) 



= -V ■ ^(1) + (V • - V ■ _ i (V . ^(i))2 - d.^'^j^^d^^ 



V-*(2) + C'(3) 



(1) 



0(3) 



(A3) 



where we have expanded in the same way that we expanded 0, i5 and v in Section |TT1 With equations (jA2p and 
[3]) . we can find the same physical quantities in Lagrangian perturbation theory (LPT) as we calculated directly in 
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Eulerian perturbation theory (EPT) in Section [III The advantage is that, in LPT, the transition to redshift space is 
merely a coordinate transformation. In this section, we show that our EPT results are consistent with LPT results 
to second order. In reference [s^l, Bernardeau claims that this is not the case, that LPT and EPT do not agree and 
that LPT is more consistent with the results of N-body simulations (see the discussion following his Figure 48). This 
is merely a confusion of terminology. What Bernardeau meant to claim is that second order LPT and EPT results 
differ substantially from third order LPT results (the difference between "tree-level" and "one-loop level" perturbation 
theory in reference [46}) and that N-body simulations are more consistent with the third order results. In Section [TVl 
we found that the effect of wq on the second order results is not strong enough to justify extending our calculations 
out to higher order. 

Before we begin, there is a subtlety which bears discussion. All of the ^ functions and their spatial derivatives in 
equations (|A2[) and (|A3p are expressed in terms of q coordinates, that is, in terms of the spatial positions of fluid 
elements at initial time r = 0. When wc go to derive field equations for and ^, we will be expressing and its 
spatial derivatives in terms of x coordinates (spatial positions now). Assume that <^ qj. Because we will only 
be interested in quantities up to second order, the transformation between q and x coordinates is a straightforward 
Taylor expansion. If we use ^!{x) to represent the displacement field as a function of the present position x (see 
reference [41|), then 

^{q) = l(f)+5,J(.f)((7-f), + 0(|3) 

= 1(f) - ^A^^ix) + e»(^3) (A4) 
dg^i'l) = |^a,(l(x)-^'AJ(f) + 0(|3)) 

where we have used the zeroth order equation q~x. Now 

5 = -V • + (I • V)V ■ I + (V ■ *(i))2 - i((V • *(i))2 - dM^dj^'i^^) - V • + 0(3) (A6) 

= - ^-f^S.^f^' -f 0(3). (A7) 

From here on, we will deal only in ^!{x) and derivatives with respect to x. Let us now consider the LPT equations of 
motion. 

LPT proceeds from the geodesic equation 

U°'VaU' = 

{i is a spatial index) and the Poisson equation. Wc once again replace the Poisson equation with the space-time 
Einstein equation in the form (jl3p 



V" (^0 + ^(1 + w)<j)j = -n'n^s 

with 5 taken from equation (|A6p so that, to first order, 

^7t:2an(v-|(i)) =-W^(^^ + n{l + w)c^). (A8) 
Using metric ([7]) and four- velocity ([S]), the geodesic equation taken out to zeroth order in 1/c^ merely returns equation 

m 

^+'Hv + {v-V)v + V(t){l+w) = 0. 
To first order in {(/>, 5, this equation gives 
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where we have used equation (|A7|) . The time derivative of equation (jA8p gives 

3, 



2 



3 



H21]™(-V20(1+ti7))+2HV2 (0 + ^(1 + ^7)0) = (AlO) 



where, in the last hne, we substitute from equations (|A8[) and (|A9p . Equation (|A10p can be rearranged to give the 
first order EPT equation (fT5|) . Clearly, the first order result for (j> is identical in either LPT or EFT formalism. The 
same is true of S and v, given that equations (jA6p and (jA7p give (5 = — V • t;, which is the first order part of the EPT 
equation of motion Q. To second order, we illustrate the consistency between LPT and EPT by observing that the 
second order parts of equations (jA6P and (jA7P are such that 

= -a,i)f' - (^(I'a.wf ' - (All) 

which is the second order part of equation ^ . Note that we have used the assumption that the fluid flow is irrotational 
so that di'^j = dj^i. Using result (jAlip as motivation, we will now write the growth functions from reference [4l| in 
terms of our own and show that our second order solutions are, indeed, equivalent to those of LPT. 

Results in reference [4l| are presented in terms of the growth factors gi and 92 and their dimensionless derivatives 
fi = {a / gi){dgi / da) . These growth factors are defined such that 

^ = .gil-f^' + 52 

where ^ denotes the spatial part of ^ . From the first order part of equation (|A3p . it is easy to see that gi is just our 
growth factor D from equation p9|) . To find (72 as a function of our second order growth functions P° and , we 
compare our equation (PT|) to equation (A19) of reference [4l[, which finds that 5 goes as 

= j '^^^[ni'hlklwi\k, + k.\y^^^^^^^ (A12) 



The comparison with our equation (|4ip proceeds more directly if we write 

2?" = -2Dl + ^0-2 

where Di is just the first order growth function D from equation p9p . and D2 is defined in equation (39a) of reference 
[47! (this equivalence can be seen by comparing that work's equation 48 with our equation I4ip . Now, we can rewrite 
our equation as 

X {{-2Dl + Id2){1 + ^ cos 012) + {^Dl - -D2){1 + cos^ 6,2 + 2^ cos0i2)) . (A13) 
\ 2 fc2 2 4 /C2 / 

Comparing equations (|A12[) and (|A13|) . we find that Di = gi (as promised) and 

P2 = 2Dl - ^D2 

or 

32 = -3P" - 41?^ + 2DI (A14) 
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0.002 




0.2 0.4 0.6 0.8 1 

e/7t 

FIG. 5: We plot the relative disagreement in Bavg between Eulerian and Lagrangian (reference [4ll. l48l|) perturbation theories. 
The horizontal axes are the angle of separation between the vectors ki and k2 in units of tt. The vertical axes are (Bavg, ept — 
Savg, LPT)/Bavg,LPT. We Compare (|55|l with 6i = 1 and 62 = to the results presented in equations (33)-(35) of reference [4l[. 
We use equation HA14|I to express g2 and /2 = {a/ g2){dg2/da) in terms of o ur g rowth functions. For all curves, zuo ~ and the 
background cosmology is the WMAP 5-year maximum likelihood universe [3a Is^. 



(A15) 
(A16) 



Using these results and equation 



we find for EPT 



V • V- 



(2) 
EPT 



(.92 + DD){V^ipV^(p + diipdiV^ip) - (32/2 + i:'I))(VVvV + di-jipdijip + 2djipdj\'^ip) 



92 



For LPT, we note that S'^^^ 
that 



^ means that * 



l^-* = —Ddiif. This, combined with the resuh from reference [4l| 

(A17) 



gives 



^ 9i 9i 



■ 92 



From this we see that equation (jA14[) . which guarantees (Slpt ~ '^ept ^^^'^ provides i^pT = ^pt- Figure [5] compares 



the results of Section [TVI with the results of LPT in the case of = 0. LPT results arc evaluated from Hivon et al's 
equations (33)-(35) using our equation (|A14p . Wc find that our EPT results agree with Hivon et al's LPT results at 
the sub-percent level. 

All of the discussion in this Appendix has assumed zu — 0. As shown above, our results arc perfectly consistent with 
LPT in this limit. One can move between the two formalisms using the algebraic relationships (|A15p and (|A16p . This 
is not the case once wq ^ 0. Equations (jAlsp and (|A16[) imply a relationship between P", V^, and D. Specifically, 
they imply 



(A18) 



Figure [H] plots the departure from this relationship as a function of for cosmologies with different values of ^Im- 
All of the curves are evaluated at redshift z = 0. While wc see that equation (|A18[) is identically true for wq = 0, 
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redshift: z = 

0.06, 1 ^ 1 ^ 1 r 



0.05 - 




2 4 6 



^0 

FIG. 6: We plot the departure of C from the predictions of equation (|A18|I as a function of vjq for different fim cosmologies. 
While equation (|A18P holds exactly when wo = 0, the relationship breaks down as \zuo\ grows. 



it becomes less and less true as |n7o| grows. Recall that equations (|A15[) and (|A16p were derived by assuming that 
^LPT — ^EPT ^^'^ comparing the growth functions of the different spatially dependent terms. If we assume that LPT 
and EPT remain in agreement for wq ^ 0, it seems likely that the spatial dependence of S^"^^ must change for vjq ^ 0. 
Specifically the spatial relationship (|A17p may break down, spoiling the relationship between our EPT growth factors 
and the LPT growth factor 172 • Given our success at reproducing the zuq = solutions of LPT, the main body of this 
paper analyzed the wq 7^ bispectrum using our EPT solutions. It may be worth solving the exact wq ^ LPT 
equations of motion in some future work. 



APPENDIX B: THE SKEWNESS 



The third order moment (the skewness) is a measure of the asymmetry of the S distribution function once non-linear 
growth has set in. Mathematically, it is defined as 



(<52(f))2 



(Bl) 



Kamionkowski and Buchalter derive an exact expression for the skewness of the matter distribution in a general 
relativistic universe under different assumptions of flatness and acceleration or deceleration [3^ . In this section, we 
follow their lead and derive a similar expression in the case of 1370 7^ 0. 

As with the bispectrum, the numerator of the skewness (jBl[) reduces to 3{5^^^{S'^^^)'^). Equation (|4ip gives us 



(J(2)(^(l))2^ 



d^ki...(Pk 



4glS-(fei+fc2+fc3+fc4) 



D^kiki 



ks ■ k^kiip" + 2©") + kikiiV + V") + {k^ ■ kifv'' 



Using the following results 

d^k,d^k2Ps{h)Ps{k2) I ,,2 



(27r)6 



1)2 



.(-2fc2 



^1-^2 — 1 



d^kid^k2 Ps{ki)Ps{k2) , ^(fci ■ fca) 



(2^) 



D 



d^k,d^k2 Ps{k^)Ps{k2) 

(27r)6 £)2 



d(cos 6) cos2 



1 + 2 COS^ 012 



d^kid^k2Psik,)Ps{k2) 

(27r)6 1)2 

f d^k,d^k2Ps{kl)Ps{k2) , ,_2a 

J — — 

2 
3 

5 f d2fcirf3fc2 Ps{k,)Ps{k2) 

3 J (27r)6 



(B2) 



(B3) 



(B4) 
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FIG. 7: We plot the change in the skewness S3 (|B1|) resulting from varying f2m and ujo- When flm is varied, too = 0. When 
vjo is varied, f2m = 0.256. AS3 is assessed relative to a vjq = 0, Qm = 0.256 universe. All other parameters are set to be 
the WMAP 5-year maximum likelihood values [sj- Over regions of interest, the effect of varying -ajo is much stronger than 
the effect of varying Q.m- Unfortunately, this will not hold when the calculation is altered to take into account a realistic 
observational window function (see Figure |8}. 



we find that the normalized skewness may be written 



((J(l))2}2 1)2 ^ 3 

Note that the denominator of the left hand side of equation (jB5P is effectively 



(27r)e 



-Ps{ki)Psih 



2) 



hence the simplicity of the right hand side. 

Figures 7(a) and |7(b)] plot the effect of the parameters flm and wq on the skewness as determined by equation (|B5 



When Q,n is varied, wq = 0. When vjq is varied, = 0.256. All other parameters are set to be the WMAP 5-year 
maximum likelihood values (33| . As in Figure [3l the vertical axis is 

AS's 5*3 — S'3|roo=0.0,„=0.256 

a = ^ ■ (B6) 

"JS.GR 'J3|roo=0,f2™=0.256 

Varying vjq within the WMAP 5-year 2a limit {wq = 1-4^2 q (2c) [l^) results in a few percent variation in 5*3, while, 
for Wq = 0, variations of flm within the WMAP 5-year 2a range result in only a few tenths of a percent change in 
S3. The la constraint reported by the WMAP 5-year team is = 0.26 ± 0.03 [3^. Though a few percent is still far 
too fine a variation to expect to detect in the data, the fact that the change due to wq is so much more pronounced 
than the change due to flm leaves open the hope that we may one day be able to distinguish between an alternative 
gravity and a dark-energy dominated universe on the basis of the matter overdensity distribution. Unfortunately, as 
we will see below, even this hope evaporates in the face of considerations associated with the mechanics of making an 
observation. 



1. The smoothed skewness 



As with the bispectrum, we cannot achieve perfect knowledge of the skewness at any scale. The birdshot distribution 
of our mass gauges (galaxy dynamics, lensing of background galaxies, etc.) forces us to accept that our measurement 
of 5{x) will have to be averaged over some scale Rq (e.g. as is the rms average mass fluctuation within an 8/i^^Mpc 
sphere). Bernardeau derives expressions for both the skewness and the kurtosis smoothed over some scale Rq in 
reference In this section, we will follow his calculation, convolving our results from Section [Bl with a top- hat 

window function. We find that whatever hope Figures [7(a)] and [7(b)] give us for constraining wq vanishes behind the 
limitations of actual measurement. 
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For the purposes of this section, we will work with a simple top-hat window function ORa{r) which is zero ii r > Rq 
and unity otherwise. The smoothed first order overdensity is therefore 



W{kRo) = 3 



d^x'd^k 
(27r)3/2 
d^x'd^k 
(27r)3/2 

d^k 
(27r)3/2 

sin(fci?o) cos(fci?o) 



(B7) 



e'^-^'~5'^^\k)QR„{\x' ~x\) 
S'^^\k)W{kRo)e'^-^ 



(B8) 



{kRo 



{kRoY 



W{kRo) is the Fourier transform of Qb„. The second order smoothed overdensity is 



d^kid^k2 ^(^)(fci)^(^)(fc2 

(27r)3 ZP 



^M-(/ci+fe2) 



W^(i?0|fcl +/C2 



fc2 



kj 



(B9) 



which is what we would expect from equation ([1T|) written in a slightly different form and subjected to the same 
convolution explicitly shown going from equation jBTj to equation (|B8p . 
By analogy with equation (jBip . the smoothed skewness is defined as 



'S'3(-Ro) 



(BIO) 



From equations (jBSp and (jB9p . we can write (note that, as in Bernardeau, Wi = W{kiRo) and Wi+j = W(|fci + £, |i?o)) 

(i^fc;^ . . . d^ki 



(27r)6 ''''l''''2^'l'3+4^2 



I pa 



^3 ■ ^4 



1 + 2 



ka -ki {ks ■ kif 



P/r2 

3 4 



((5i(52(53(54)- 



(Bll) 



Remembering our power spectrum convention (|40)) . it is useful to note that the coefficients of both the 2?° and the 
T)^ terms will vanish if fc3 ~ —k^. Therefore, the (5i(52)(^3'54} term does not contribute to the integral and we have 



= 6 



d^ki...d^k. 



i2n) 



g WiW2W3+i^ 



fc3 • h 



V" I 1 + 2 



^3-^4 , (^3 • h)'^ 



i.2i.2 



X Ps{ki)Ps{k2)SUki + fc3)4(fc2 + ki) 



d^kid^k2 1 
^^W,W2W,^2^^ 



xP(fci)P(fc2) 



ki ■ k2 



1 + 2 



ki ■ k2 {ki ■ fc2)^ 



fc2 



k^k^ 
"-1 "'2 



(B12) 
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FIG. 8: We plot the smoothed skewness (|B10P as a function of scale Ro for different values of roo and Qm- Results are presented 
(as in Figure [T} as a change in skewness relative to the WMAP 5-year maximum likelihood GR cosmology. While zuo retains 
the strong influence on S3 observed in Figure 7(b) the effect of Qrn on the shape of the power spectrum amplifies its infiuence, 
so thtat it is now very difficult to distinguish between a significant departure from general relativity and a minor shift in Qrn- 
The la confidence interval reported by the WMAP 5-year team is fim ~ 0.26 ± 0.03. 



It is tempting to proceed as in Section [B] and analytically evaluate equation (|B12p to find the same result as equation 
(|B5I) with an added factor of 

(ZTT j" 

Unfortunately, the compound window function W1+2 depends on the angle between ki and fc2, and thus spoils the 
straightforward angular integrations leading to equations ()B3| and (|B4p . Furthermore, the denominator of equation 
(|B10[) has a completely different dependence on the window functions and power spectra. To wit 



2\2 



(27r)3 



Ps{k)W^ikRo) 



(B13) 



The smoothed skewness (jBlOp cannot be written as a simple function of growth factors. It must be integrated over 
the 6 dimensional parameter space {ki, fc2}, which can be reduced to a 4 dimensional parameter space upon noting 
that the integrand in equation (|B12p does not explicitly depend on the azimuthal angles (f>ki or so that we can 
replace those two integrations with a factor of 47r^. We use the same Monte Carlo integrator wc used in Section [TVl 
We present our results in Figure (jS]) as a change in skewness (relative to a WMAP 5-year maximum likelihood GR 

the effect of flm 



universe; see equation IB6|) . While vdq retains the strong influence on 6*3 observed in Figure 7(b) 
on the shape of the power spectrum amplifies its influence, so that it is now very difficult to distinguish between a 
significant departure from general relativity and a minor shift in Qm- Wc perform a similar calculation for 6*4 (i?o) in 
Section O 



APPENDIX C: FOURTH ORDER MOMENTS 



The smoothed kurtosis 54(i?o) of the overdensity distribution is the fourth order moment of the distribution averaged 
over a radius Rq, i.e. 



S4{Rf)) — 



(CI) 



In this section, we will show the effect of 7^ on 6*4 (i?o)- While the effect is stronger than that found on the 
skewness in Section |Bl this is of little comfort because, as a higher order moment, the kurtosis is much harder to 
measure than the (already difficult) skewness. Also, we find that the effect of varying fl„i is similarly amplified so 
that it would still be difficult to discern a large change in zuq from a small change in Qm- 
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1. Third order equations of motion 



Calculating the kurtosis requires first that we find the evolution of the third order perturbations {if>^^\S^^\v''^'>}. 
We do so in this section. For easy reference, we recall the first and second order solutions 





= /(r)(^(f) 


^(1) 


= D{t)V^lp 




= -DVip{x) 




= a{T)A{x) + P{t)B{x) 


^(2) 


= P(")(t)V2A + (r)V2B 




= + Db)VA - V^VB 











The third-order part of equation (|12p is 

j(3) + ^^(3) ^ (1 ^ ^) (a,<5(2)5^^(i) ^ a^^(i)a^0(2)) + g^^. [2i;*(2)^,(i) ^ ^(i)„.(i)^,(i) 

+ (1 + n7)(<5(i)v20(2) + 5(2)v20(i) + V20(3)). 
Using the recipe outlined in Section [TTl we find 

0(3) ^ ^F(')(r)e«(a?) 



(C2) 



7 



J(3) 



(C3) 



where 



v2g(«) = d,,[d,^d,B] 

The equations of motion for i^^*'' and J^*-*^ are 

pi^) = -F(')7i(3 + Tu) - [^^3(1 + n7)(l - f7„0 
i^'') = -WJ^''') + (l + n7)F('''+S'(3'''). 
The right hand side of equation (|C2p gives 

5^3'^' = (l + n7)/2?('^) 
5(3'^' = {l + m)Da 
^(3,5) ^ -2b{-V'' + DD) 



(C4) 



5(3,6) 
5(3.7) 
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2. The kurtosis 



The numerator of equation (jCl[) is equivalent to 
Using equation (|C3p and drawing an analogy with equation (jB12[) . we find that 



^,a)3,(3)^ ^ 4 f d^k^dWs P(fcr)P(fc.)P(fc3) ^^^^^3^^^^^3 ^ ^(.) ^ ^«(fc,fc,,,^3) (C5) 

^ perm 



where X^pcrm^^*^ = ^2, ^3) + ^^'H^ii ^3^ ^2) + ••■ and the are the Fourier transforms of equations (|C4p 

with all factors of 1^9 artificially removed. 

Similarly, the other part of the kurtosis reduces to 

^,(r)2,(2)2^ ^ , J d3W3fe3p(fe^)p(,^)p(,3) ^^^^ ^ W^3-rM..3-.l3(fc3, -k^M^k,, ^k.) (C6) 

(C7) 



perm 



In this case X^porm W^3- 1^^-3-22^(^3, -~ki)V{—k3, —k-i) has 8 terms, permuting over both the content of the two pairs 
of vectors 

P(fc3, -ki)V{-k3, -k2) V{k3, -fc2)P(-A?3, -fci), etc. 
and permuting over the vectors' orders within the pairs 

V{k3, -ki)T){-k3, -k2) V{-ki,k3)f){-k2, -fcs), etc. 

Note also that fca is always positive in the first pair and negative in the second, fci and ^2 are always negative. These 
minus signs are acquired from the delta function in 

{5^^\h)~5^^\kj))=P{h)5l{h + kj). 

Figure [H recreates Figure [5] for the smoothed kurtosis. Again wc find that it is nearly impossible to distinguish a large 
change in voq from a small change in 17^- 
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FIG. 9: We plot the smoothed kurtosis HC1|I as a function of scale Ro for different values of wo and Q,m- Results are presented 
(as in Figure [3| as a change in kurtosis relative to the WMAP 5-year maximum likelihood GR cosmology. The la confidence 
interval reported by the WMAP 5-year team is Q,m = 0.26 ± 0.03. 



